#####################
### Load Packages ###
#####################

library(haven)


# Set this Working Directory to the Folder with the Replication Files
#setwd("C:/Steven/GDrive/Research Files/Clark-Rogers/Descriptive Representation/APSR Final Submission/Replication")

###
# Functions
###

# Confidence Interval for Weighted Means
# https://stackoverflow.com/questions/37973240/how-to-create-a-confidence-interval-for-a-weighted-average-of-areas-under-the-ro
weighted.ttest.ci <- function(x, weights, conf.level = 0.95, na.rm=TRUE) {
  #x <- SLUCCES$Know_TF_Abortion_Correct
  #weights <- SLUCCES$teamweight
  require(Hmisc)
  nx <- length(x[!is.na(x)]) # Rogers editted here to account for missing data
  df <- nx - 1
  vx <- wtd.var(x, weights, normwt = TRUE) ## From Hmisc
  mx <- weighted.mean(x, weights, na.rm=TRUE)
  stderr <- sqrt(vx/nx)
  tstat <- mx/stderr ## not mx - mu
  alpha <- 1 - conf.level
  cint <- qt(1 - alpha/2, df)
  cint <- tstat + c(-cint, cint)
  cint * stderr
  
}

# Create Dataframe of CIs
CIDataframe <- function(CurrentVariable, CurrentWeight, Issue)
{
  WeightedCI <- weighted.ttest.ci(CurrentVariable, CurrentWeight)
  CurrentWeightedCIs <- data.frame(matrix(NA, 1,0))
  CurrentWeightedCIs$Issue <- NA
  CurrentWeightedCIs$LowerBound <- NA
  CurrentWeightedCIs$UpperBound <- NA
  
  CurrentWeightedCIs$Issue[1] <- Issue
  CurrentWeightedCIs$LowerBound[1] <- WeightedCI[1]
  CurrentWeightedCIs$UpperBound[1] <- WeightedCI[2]
  return(CurrentWeightedCIs)
}

ClarkRogersData <- read_dta("Data/Formatted Data/Clark-Rogers-RegressionData.dta")

# Create Subsets

AllVoters <- ClarkRogersData
FemaleVoters <- subset(ClarkRogersData, ClarkRogersData$Female==1)
BlackVoters <- subset(ClarkRogersData, ClarkRogersData$Black==1)
CollegeDegree <- subset(ClarkRogersData, ClarkRogersData$CollegeDegree==1)

BarPlotBin <- matrix(NA,2,0)
  BarPlotBin$WomenCongress <- NA
  BarPlotBin$WomenCongress[1] <- weighted.mean(ClarkRogersData$WomenCongress_Correct,ClarkRogersData$teamweight, na.rm=TRUE)
  BarPlotBin$WomenCongress[2] <- weighted.mean(FemaleVoters$WomenCongress_Correct,FemaleVoters$teamweight, na.rm=TRUE)
  BarPlotBin$WomenCongress[3] <- weighted.mean(BlackVoters$WomenCongress_Correct,BlackVoters$teamweight, na.rm=TRUE)
  BarPlotBin$WomenCongress[4] <- weighted.mean(CollegeDegree$WomenCongress_Correct,CollegeDegree$teamweight, na.rm=TRUE)

  BarPlotBin$WomenStateLeg <- NA
  BarPlotBin$WomenStateLeg[1] <- weighted.mean(ClarkRogersData$WomenStateLeg_Correct,ClarkRogersData$teamweight, na.rm=TRUE)
  BarPlotBin$WomenStateLeg[2] <- weighted.mean(FemaleVoters$WomenStateLeg_Correct,FemaleVoters$teamweight, na.rm=TRUE)
  BarPlotBin$WomenStateLeg[3] <- weighted.mean(BlackVoters$WomenStateLeg_Correct,BlackVoters$teamweight, na.rm=TRUE)
  BarPlotBin$WomenStateLeg[4] <- weighted.mean(CollegeDegree$WomenStateLeg_Correct,CollegeDegree$teamweight, na.rm=TRUE)

  BarPlotBin$BlackCongress <- NA
  BarPlotBin$BlackCongress[1] <- weighted.mean(ClarkRogersData$BlackCongress_Correct,ClarkRogersData$teamweight, na.rm=TRUE)
  BarPlotBin$BlackCongress[2] <- weighted.mean(FemaleVoters$BlackCongress_Correct,FemaleVoters$teamweight, na.rm=TRUE)
  BarPlotBin$BlackCongress[3] <- weighted.mean(BlackVoters$BlackCongress_Correct,BlackVoters$teamweight, na.rm=TRUE)
  BarPlotBin$BlackCongress[4] <- weighted.mean(CollegeDegree$BlackCongress_Correct,CollegeDegree$teamweight, na.rm=TRUE)
 
  BarPlotBin$BlackStateLeg <- NA
  BarPlotBin$BlackStateLeg[1] <- weighted.mean(ClarkRogersData$BlackStateLeg_Correct,ClarkRogersData$teamweight, na.rm=TRUE)
  BarPlotBin$BlackStateLeg[2] <- weighted.mean(FemaleVoters$BlackStateLeg_Correct,FemaleVoters$teamweight, na.rm=TRUE)
  BarPlotBin$BlackStateLeg[3] <- weighted.mean(BlackVoters$BlackStateLeg_Correct,BlackVoters$teamweight, na.rm=TRUE)
  BarPlotBin$BlackStateLeg[4] <- weighted.mean(CollegeDegree$BlackStateLeg_Correct,CollegeDegree$teamweight, na.rm=TRUE)   
 
  BarPlot <- cbind(BarPlotBin$WomenCongress, BarPlotBin$WomenStateLeg, BarPlotBin$BlackCongress, BarPlotBin$BlackStateLeg)
  
  
    ###
    # Confidence Intervals
    ###
    
      CIs <- CIDataframe(ClarkRogersData$WomenCongress_Correct,ClarkRogersData$teamweight, "Women-Congress-All")
      CIs <- rbind(CIs, CIDataframe(FemaleVoters$WomenCongress_Correct,FemaleVoters$teamweight, "Women-StateLeg-All"))
      CIs <- rbind(CIs, CIDataframe(BlackVoters$WomenCongress_Correct,BlackVoters$teamweight, "Black-Congress-All"))
      CIs <- rbind(CIs, CIDataframe(CollegeDegree$WomenCongress_Correct,CollegeDegree$teamweight, "Black-StateLeg-All"))
      
      CIs <- rbind(CIs, CIDataframe(ClarkRogersData$WomenStateLeg_Correct,ClarkRogersData$teamweight, "Women-Congress-Female"))
      CIs <- rbind(CIs, CIDataframe(FemaleVoters$WomenStateLeg_Correct,FemaleVoters$teamweight, "Women-StateLeg-All"))
      CIs <- rbind(CIs, CIDataframe(BlackVoters$WomenStateLeg_Correct,BlackVoters$teamweight, "Black-Congress-All"))
      CIs <- rbind(CIs, CIDataframe(CollegeDegree$WomenStateLeg_Correct,CollegeDegree$teamweight, "Black-StateLeg-All"))
      
      CIs <- rbind(CIs, CIDataframe(ClarkRogersData$BlackCongress_Correct,ClarkRogersData$teamweight, "Women-Congress-Female"))
      CIs <- rbind(CIs, CIDataframe(FemaleVoters$BlackCongress_Correct,FemaleVoters$teamweight, "Women-StateLeg-All"))
      CIs <- rbind(CIs, CIDataframe(BlackVoters$BlackCongress_Correct,BlackVoters$teamweight, "Black-Congress-All"))
      CIs <- rbind(CIs, CIDataframe(CollegeDegree$BlackCongress_Correct,CollegeDegree$teamweight, "Black-StateLeg-All"))

      CIs <- rbind(CIs, CIDataframe(ClarkRogersData$BlackStateLeg_Correct,ClarkRogersData$teamweight, "Women-Congress-Female"))
      CIs <- rbind(CIs, CIDataframe(FemaleVoters$BlackStateLeg_Correct,FemaleVoters$teamweight, "Women-StateLeg-All"))
      CIs <- rbind(CIs, CIDataframe(BlackVoters$BlackStateLeg_Correct,BlackVoters$teamweight, "Black-Congress-All"))
      CIs <- rbind(CIs, CIDataframe(CollegeDegree$BlackStateLeg_Correct,CollegeDegree$teamweight, "Black-StateLeg-All"))
        
    
    # Find Bar Centers
    
      
      
    BarCenters <- barplot(BarPlot,
                          names.arg=c(" ", " ", " ", " "), beside=TRUE,
                          ylim=c(0, 1),
                          yaxt="n",
                          cex.names =.75,
                          ann=FALSE
                          )
    
                  BarCenters <- t(as.data.frame(BarCenters))
                  CIs$BarCenters <- as.numeric(c(BarCenters[1,], BarCenters[2,], BarCenters[3,], BarCenters[4,]))
                  dev.off()
    
    ###
    # Plot Time
    ###
    
    tiff(filename = "Figures/Figure-3.tiff", width = 7, height = 5, units = "in", pointsize = 12, res=1000)
                  
                  
                  windowsFonts(
                    A=windowsFont("Arial Black"),
                    B=windowsFont("Cambria"),
                    C=windowsFont("Comic Sans MS"),
                    D=windowsFont("Symbol")
                  )
     
    par(family="B",
        mai=c(.5,1,.25,.1))
                               
    barplot(BarPlot,
            names.arg=c(" ", " ", " ", " "), beside=TRUE,
            ylim=c(0, .55),
            yaxt="n",
            cex.names =.75,
            ann=FALSE
            )
      
    axis(2, at=c(.0,.1,.2,.3,.4,.5),
         labels=c("0%", "10%", "20%", "30%", "40%", "50%"), cex.axis=1, las=2)
      
    axis(1, at=c(3, 8, 13, 18),
         labels=c("Women", "Women", "Black", "Black"), las=0, tick=FALSE, line=-1, cex.axis=.9)
    axis(1, at=c(3, 8, 13, 18),
         labels=c("Congress", "State Leg.", "Congress", "State Leg."), las=0, tick=FALSE, line=0, cex.axis=.9)
    
    
    par(cex.axis=1.5)
    
    axis(2, at=c(.25),
         labels=c("Percentage Correct Responses"), cex=6, line=2.2, cex.axis=1.3, tick=FALSE)
    
    colors <- gray.colors(4) # Number of Gray Colors Used
    
    # Add Confidence Intervals
    segments(CIs$BarCenters, CIs$LowerBound, CIs$BarCenters, CIs$UpperBound, lwd = 1.5)
    arrows(CIs$BarCenters, CIs$LowerBound, CIs$BarCenters, CIs$UpperBound, lwd = 1.5, angle = 90, code = 3, length = 0.05)
    
    
    legend(5,.55,
           c("All Respondents", "Women Respondents", "Black Respondents", "College Degree Resp."),
           fill=colors,
           ncol=2,
           cex=1)
    
    dev.off()
    